

** make bootstrap cis
/*
1. append together bootstrap sample
2. in each sample-week, calcualte bounds based on non-icli and clear cause hospitalizations
3. find 2.5th and 97.5th percentile of bounds (by week) 
4. find average lb and ub, by week, and get weekly bias in bound 
*/

** 1. append together
local B 500
clear
forvalues bs = 1/`B' {
	append using bootstrap_ests/bs_u_`bs'
}

** 2. calculate intersection bounds 
gen lb_non_m = lower_bound1
gen lb_clear_m = lower_bound1
gen lb_icli_m = lower_bound1 


gen lb_non_i = max( lower_bound1, lower_bound2)
gen lb_clear_i = max(lower_bound1, lower_bound3)
gen lb_icli_i = max(lower_bound1 , lower_bound4)

gen ub_non_m = min(upper_bound1, upper_bound2)
gen ub_clear_m = min(upper_bound1, upper_bound3)
gen ub_icli_m = min(upper_bound1, upper_bound4)

gen ub_non_i = ub_non_m 
gen ub_clear_i = ub_clear_m  
gen ub_icli_i = ub_icli_m 

** 3. find percentiles
local p025 = ceil(`B'*.025)
local p975 = ceil(`B'*.975)
foreach t in non clear icli {
	foreach suff in i m {
		sort week lb_`t'_`suff'
		by week: gen p025_`t'_`suff' = lb_`t'_`suff'[`p025']
		sort week ub_`t'_`suff'
		by week: gen p975_`t'_`suff' = ub_`t'_`suff'[`p975']
	}	
}
sort week lower_bound1
by week: gen p025_pop = lower_bound1[`p025']
sort week upper_bound1
by week: gen p975_pop = upper_bound1[`p975']

** find means
foreach t in non clear icli {
	gen flip_`t' = ub_`t'_m == upper_bound1
}
sum flip_* 

sort week 
collapse (mean) lb_non_m-ub_icli_i flip* p025_* p975_*, by(week) 


foreach t in non clear icli {
	foreach b in lb ub {
		foreach suff in i m {
			rename `b'_`t'_`suff' mean_`b'_`t'_`suff'
		}
	}
}
save tables/bootstrap_ci_unwtd, replace

** set up point estimates and cis   estimates 
use bootstrap_ests/bs_u_0 , clear

gen lb_non_m = lower_bound1
gen lb_clear_m = lower_bound1
gen lb_icli_m = lower_bound1 


gen lb_non_i = max( lower_bound1, lower_bound2)
gen lb_clear_i = max(lower_bound1, lower_bound3)
gen lb_icli_i = max(lower_bound1 , lower_bound4)

gen ub_non_m = min(upper_bound1, upper_bound2)
gen ub_clear_m = min(upper_bound1, upper_bound3)
gen ub_icli_m = min(upper_bound1, upper_bound4)

gen ub_non_i = ub_non_m 
gen ub_clear_i = ub_clear_m  
gen ub_icli_i = ub_icli_m 

rename lower_bound1 lb_pop
rename upper_bound1 ub_pop

merge 1:1 week using tables/bootstrap_ci_unwtd, assert(3) 
drop _merge

** bias ub m 
foreach t in non clear icli {
	gen bias_`t' = ub_`t'_m - mean_ub_`t'_m
}

save tables/estimates_unwtd_bs_ci, replace 

